
library(haven)
library(devtools)
install_github("pewresearch/pewmethods", build_vignettes = TRUE)
library(pewmethods)
library(tidyverse)
dat <- read.csv("~/Desktop/jop replication/sample_for_weighting.csv")

# These need to be specified for the User
outdir_tables <-"~/Desktop/jop replication/"
outdir_figures <-"~/Desktop/jop replication/"


# 1 = Donor, 2= High Income, 3 = Gen Pop

##########################################################################################
#     All Donors
###########################################################################################

dat$resp.indx <- seq(1,nrow(dat))
rownames(dat) <- dat$resp.indx

dat.donor <- dat[dat$sample==1,]

dat.donor$ts_genderstring[is.na(dat.donor$ts_genderstring)] <- ""
dat.donor$ts_race_black[is.na(dat.donor$ts_race_black)] <- ""
dat.donor$tsd_registeredvote[is.na(dat.donor$tsd_registeredvote)] <- ""
dat.donor$tsd_regdem[is.na(dat.donor$tsd_regdem)] <- ""
dat.donor$tsd_regrep[is.na(dat.donor$tsd_regrep)] <- ""
dat.donor <- dat.donor[-as.logical(dat.donor$tsd_regdem=="")*seq(1,nrow(dat.donor)),]

nrow(dat.donor)
sum(dat.donor$completedsurvey,na.rm=TRUE)
mean(dat.donor$completedsurvey[dat.donor$tsd_regrep=="yes"],na.rm=TRUE)
mean(dat.donor$completedsurvey[dat.donor$tsd_regdem=="yes"],na.rm=TRUE)

age.quartiles <- summary(dat.donor$ts_age)
age.target <- age.quartiles[c(2,3,5)]

dat.donor$w.age <- NULL
dat.donor$w.age[dat.donor$ts_age < age.target[1]] <- "<53"
dat.donor$w.age[dat.donor$ts_age >= age.target[1] & dat.donor$ts_age < age.target[2]] <- "53-64"
dat.donor$w.age[dat.donor$ts_age >= age.target[2] & dat.donor$ts_age < age.target[3]] <- "64-73"
dat.donor$w.age[dat.donor$ts_age >= age.target[3]] <- "73-100"
dat.donor$w.age[is.na(dat.donor$ts_age)] <- "Missing"
dat.donor$w.age <- as.factor(dat.donor$w.age)

# Do for partisan score

partisan.target <- summary(dat.donor$ts_tsmart_partisan_score)[c(2,3,5)]
dat.donor$party.score <- NULL
dat.donor$party.score[dat.donor$ts_tsmart_partisan_score <= partisan.target[1]] <- "<=4"
dat.donor$party.score[dat.donor$ts_tsmart_partisan_score > partisan.target[1] & dat.donor$ts_tsmart_partisan_score < partisan.target[2]] <- "5-66"
dat.donor$party.score[dat.donor$ts_tsmart_partisan_score >= partisan.target[2] & dat.donor$ts_tsmart_partisan_score < partisan.target[3]] <- "67-97"
dat.donor$party.score[dat.donor$ts_tsmart_partisan_score >= partisan.target[3]] <- "98+"
dat.donor$party.score[is.na(dat.donor$ts_tsmart_partisan_score)] <- "Missing"
dat.donor$party.score <- as.factor(dat.donor$party.score)

dat.donor$ts_genderstring <- as.factor(dat.donor$ts_genderstring)

dat.donor$w.ts.wealth <- NULL
dat.donor$w.ts.wealth[dat.donor$ts_wealth=="Below $25,000" | dat.donor$ts_wealth=="$25,000 - $49,999" | dat.donor$ts_wealth=="$50,000 - $74,999" | dat.donor$ts_wealth=="$75,000 - $99,999"] <- " < $100,000"
dat.donor$w.ts.wealth[dat.donor$ts_wealth=="$100,000 - $149,999" | dat.donor$ts_wealth=="$150,000 - $199,999"] <- "$100,000 - $199,999"
dat.donor$w.ts.wealth[dat.donor$ts_wealth=="$200,000 - $249,999" | dat.donor$ts_wealth=="$250,000 - $499,999"] <- "$200,000 - $499,999"
dat.donor$w.ts.wealth[dat.donor$ts_wealth=="$500,000 - $749,999" | dat.donor$ts_wealth=="$750,000 - $999,999"] <- "$500,000 - $999,999"
dat.donor$w.ts.wealth[dat.donor$ts_wealth=="$1,000,000 - $2,499,999"] <- "$1,000,000 - $2,499,999"
dat.donor$w.ts.wealth[dat.donor$ts_wealth=="$2,500,000 - $4,999,999" | dat.donor$ts_wealth=="$5,000,000 - $9,999,999" | dat.donor$ts_wealth=="Above $10,000,000"] <- ">$2,500,000"
dat.donor$w.ts.wealth[dat.donor$ts_wealth=="Missing"] <- "Missing"

# Make factors for weighting

dat.donor$tsd_regrep <- as.factor(dat.donor$tsd_regrep)
dat.donor$tsd_regdem <- as.factor(dat.donor$tsd_regdem)
dat.donor$tsd_registeredvote <- as.factor(dat.donor$tsd_registeredvote)
dat.donor$ts_race_black <- as.factor(dat.donor$ts_race_black)
dat.donor$w.ts.wealth <- as.factor(dat.donor$w.ts.wealth)
dat.donor$ts_vote_2016_presprimary <- as.factor(dat.donor$ts_vote_2016_presprimary)
dat.donor$ts_vote_2016_general <- as.factor(dat.donor$ts_vote_2016_general)
dat.donor$ts_vote_2018_general <- as.factor(dat.donor$ts_vote_2018_general)

dat.donor$TotalContrib <- NULL
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont==0] <- "0"
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont==1] <- "1"
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont==2] <- "2"
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont==3] <- "3"
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont==4] <- "4"
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont > 4 & dat.donor$ts_fec_tot_cont< 10] <- "5 - 9"
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont >= 10 & dat.donor$ts_fec_tot_cont< 20] <- "10 - 19"
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont >= 20 & dat.donor$ts_fec_tot_cont< 50] <- "20 - 49"
dat.donor$TotalContrib[dat.donor$ts_fec_tot_cont >= 50] <- "50+"
dat.donor$TotalContrib <- as.factor(dat.donor$TotalContrib)


###############################################################################################################
# Appendix Table G1: Population Targets for Appendix Table
###############################################################################################################

table(dat.donor$w.age)/nrow(dat.donor)*100
table(dat.donor$tsd_regdem)/nrow(dat.donor)*100
table(dat.donor$tsd_regrep)/nrow(dat.donor)*100
table(dat.donor$party.score)/nrow(dat.donor)*100
table(dat.donor$ts_genderstring)/nrow(dat.donor)*100
table(dat.donor$ts_race_black)/nrow(dat.donor)*100
table(dat.donor$w.ts.wealth)/nrow(dat.donor)*100
table(dat.donor$ts_vote_2016_general)/nrow(dat.donor)*100
table(dat.donor$ts_vote_2016_presprimary)/nrow(dat.donor)*100
table(dat.donor$ts_vote_2018_general)/nrow(dat.donor)*100
table(dat.donor$TotalContrib)/nrow(dat.donor)*100

# Create Weighting Targets using the created factors above

targets <- create_raking_targets(
  dat.donor,
  vars = c("w.age","ts_genderstring","tsd_regrep","tsd_regdem","ts_race_black","w.ts.wealth","TotalContrib","ts_vote_2016_presprimary","ts_vote_2018_general","party.score"),
  prefix = "",
  wt = 1,
)

dat.donor.complete <- dat.donor[dat.donor$completedsurvey==1,]
wgt <- rake_survey(dat.donor.complete,pop_margins = targets)
dat.donor.complete$RakingWgt <- wgt

###############################################################################################################
# Appendix Table G1: Sample Demographics for Appendix Table
###############################################################################################################

table(dat.donor.complete$w.age)/nrow(dat.donor.complete)*100
table(dat.donor.complete$tsd_regdem)/nrow(dat.donor.complete)*100
table(dat.donor.complete$tsd_regrep)/nrow(dat.donor.complete)*100
table(dat.donor.complete$party.score)/nrow(dat.donor.complete)*100
table(dat.donor.complete$ts_genderstring)/nrow(dat.donor.complete)*100
table(dat.donor.complete$ts_race_black)/nrow(dat.donor.complete)*100
table(dat.donor.complete$w.ts.wealth)/nrow(dat.donor.complete)*100
table(dat.donor.complete$ts_vote_2016_general)/nrow(dat.donor.complete)*100
table(dat.donor.complete$ts_vote_2016_presprimary)/nrow(dat.donor.complete)*100
table(dat.donor.complete$ts_vote_2018_general)/nrow(dat.donor.complete)*100
table(dat.donor.complete$TotalContrib)/nrow(dat.donor.complete)*100

DemDenR <- density(dat.donor.complete$RakingWgt[dat.donor.complete$OwnParty=="Democratic"])
RepDenR <- density(dat.donor.complete$RakingWgt[dat.donor.complete$OwnParty=="Republican"])

# Propensity Score
prop.score <- glm(completedsurvey~tsd_regrep+tsd_regdem+w.ts.wealth+ts_race_black+ts_genderstring + w.age + ts_vote_2016_presprimary + ts_vote_2018_general + TotalContrib + party.score,
                  data=dat.donor,family=binomial(),x=TRUE)

p.score <- round(prop.score$fitted.values,digits=4)
resp.indx <- as.numeric(names(p.score))
Prop.Dat <- as.data.frame(cbind(resp.indx,p.score))
colnames(Prop.Dat) <- c("resp.indx","PropensityWgt")

dat.donor.complete <- merge(dat.donor.complete,Prop.Dat,by="resp.indx",all.x=TRUE)

# Rescale Inverse Propensity so it sums to 2000
dat.donor.complete$c1 <- 1/dat.donor.complete$PropensityWgt
dat.donor.complete$c2 <- dat.donor.complete$c1/sum(dat.donor.complete$c1,na.rm=TRUE)
dat.donor.complete$InvPropensityWgt <- dat.donor.complete$c2 * 2000

# Plot Distribution of Weights: Inverse Propensity Score Weights

DemDenP <- density(dat.donor.complete$InvPropensityWgt[dat.donor.complete$OwnParty=="Democratic"],na.rm=TRUE)
RepDenP <- density(dat.donor.complete$InvPropensityWgt[dat.donor.complete$OwnParty=="Republican"],na.rm=TRUE)


###############################################################################################################
# Appendix Figure G1: Distribution of Iterative Raking and Inverse Propensity Score weights for Donors by Self-Identified Partisanship
###############################################################################################################

pdf(file= paste(outdir_figures, "SI_FigureG1_DonorPropensityWeightDensity.pdf", sep = ""))
par(mfrow=c(1,2))
plot(DemDenR$x,DemDenR$y,type="l",col="BLUE",
     xlab="Weight",
     ylab="Density",
     main="Iterative Raking Weight")
points(RepDenR$x,RepDenR$y,type="l",col="RED")
plot(DemDenP$x,DemDenP$y,type="l",col="BLUE",
     xlab="Weight",
     ylab="Density",
     main="Inverse Propensity Score Weight")
points(RepDenP$x,RepDenP$y,type="l",col="RED")
dev.off()


